Entropy production and equilibration in Yang-Mills quantum mechanics 
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The Husimi distribution provides for a coarse grained representation of the phase space distri- 
bution of a quantum system, which may be used to track the growth of entropy of the system. 
We present a general and systematic method of solving the Husimi equation of motion for an iso- 
lated quantum system, and we construct a coarse grained Hamiltonian whose expectation value 
is exactly conserved. As an application, we numerically solve the Husimi equation of motion for 
two-dimensional Yang-Mills quantum mechanics (the x-y model) and calculate the time evolution 
of the coarse grained entropy of a highly excited state. We show that the coarse grained entropy 
saturates to a value that coincides with the microcanonical entropy corresponding to the energy of 
the system. 



I. INTRODUCTION 

Entropy production in isolated quantum systems is an 
interesting and important research problem with many 
applications. Due to the unitarity of time evolution in 
quantum mechanics, the von Neumann entropy of an iso- 
lated quantum system remains fixed. A proper definition 
of the concept of entropy growth for an isolated quan- 
tum system thus requires coarse graining which, in turn, 
must be grounded on a correspondence between quan- 
tum and classical physics. Such a correspondence can 
be constructed from one of the phase-space representa- 
tions of quantum theories found since the classical works 
of Wigner and Moyal 0, [1] . Recently, it was suggested 
by Kunihiro et al. Q that the Husimi representation of 
the density operator is suitable for describing the 

entropy production in an isolated quantum system, be- 
cause the long-term growth rate of the entropy defined 
by the Husimi distribution approaches the classical limit 
for long times. Applications of this formalism include the 
entropy production in relativistic heavy-ion collisions and 
inflationary cosmology. 

The process of entropy production in relativistic heavy- 
ion collisions has been studied extensively. The final en- 
tropy per unit rapidity produced in high-energy nuclear 
collisions at the Relativistic Heavy-Ion Collider (RHIC) 
is well known experimentally. The final entropy produced 
per unit rapidity produced in central Au+Au collisions 
at the top RHIC energy of 200 GeV per nucleon pair in 
the center-of-mass frame is 5600 ± 500 at mid-rapidity 
[?J • Theoretical studies suggest that at least half of the 
final entropy is produced during a rapid equilibration and 
thermalization period during the initial phase of the nu- 
clear collision, with a thermalization time about 1.5 fm/c 
or less [1, H| . Furthermore, it has been pointed out that 
the nuclear matter is transformed in this rapid equilibra- 
tion stage from saturated gluonic matter in a universal 
quantum state, called the color-glass condensate, into a 
thermally equilibrated quark-gluon plasma [l(j, 111] • It is 
an important theoretical challenge to construct a formal- 
ism capable of describing the entropy production during 
this equilibration and thermalization process. 

Another interesting exploration relevant to entropy 



production of quantum systems is reheating of the uni- 
verse after inflation [l2j]. The reheating process starts 
from a preheating phase [HJ], where the inflation field is 
coupled to the matter fields and it transfers energy to the 
matter fields. These matter fields then undergo further 
decays into other particles until the decay products will 
eventually reach a thermal equilibrium. Through these 
stages, the reheating process of the universe after infla- 
tion produces a gigantic amount of entropy. 

To deal with applications to such a wide range of phys- 
ical systems, it is desirable to construct a general formal- 
ism describing the coarse grained entropy production in 
an isolated quantum system from the growth of com- 
plexity of the quantum system. In this work, we ap- 
ply the formalism developed in 3] to study the coarse 
grained entropy production of a specific non-integrable 
quantum system and its approach to thermal equilib- 
rium. As an example, we choose a simple quantum sys- 
tem that possesses chaotic dynamical behavior. It is 
well-known that chaotic dynamical behavior requires that 
an isolated, conservative dynamical system must have 
at least four degrees of freedom (two position and two 
momentum variables) The two-dimensional quan- 

tum system we have chosen, often called the xy-modcl 
or two-dimensional Yang-Mills quantum mechanics, has 
well known chaotic properties [15|. We find that the 
coarse grained entropy production of this quantum sys- 
tem saturates, and we obtain a characteristic time after 
which the complexity of the system no longer increases. 

This article is structured as follows. In Section II, we 
briefly introduce the Husimi representation of the den- 
sity operator and explain how it is applied to a defi- 
nition of the coarse grained entropy of a quantum sys- 
tem, also known as the Wehrl-Husimi entropy. On the 
way, we propose an novel method to derive the coarse 
grained Hamiltonian whose expectation value serves as 
a constant of motion for time evolution of the Husimi 
distribution. In Section III, we discuss the equation of 
motion of the Husimi distribution and introduce the test- 
particle method for obtaining the numerical solutions to 
this equation. After transforming the Husimi equation 
of motion into a system of equations of motion for test 
particles, we solve these equations to obtain the Husimi 
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distribution and the Wehrl-Husimi entropy as a function 
of time in Section IV. We analyze the time dependence 
of the Wehrl-Husimi entropy and obtain a characteristic 
time scale, after which the entropy is saturated. Besides, 
we propose a method to investigate the value of the sat- 
urated Wehrl-Husimi entropy for an infinitely large test- 
particle number, which is independent of the test-particle 
approximation scheme. Furthermore, we compare the 
saturation value of the Wehrl-Husimi entropy to equilib- 
rium based definitions of the entropy of the same quan- 
tum system. The difference between the microcanonical 
and the Wehrl-Husimi entropy serves as a probe for when 
and whether the quantum system equilibrates. 



II. GENERAL FORMALISM 
A. Husimi distribution and coarse grained entropy 

The main goal of this paper is to study the entropy pro- 
duction of a quantum system as a function of time. To de- 
fine a coarse grained entropy, it is necessary to construct 
a mapping which not only creates a correspondence be- 
tween the dynamics of the quantum system and that of 
the classical system, but also ensures that the resulting 
coarse grained distribution is non-negative and thus can 
be used for the definition of the coarse grained entropy 
Q. A minimal coarse graining of a quantum system is 
achieved by projecting its density operator on a coherent 
state [4|. The resultant distribution function is known 
as the Husimi distribution pn(t;q_, p), which is positive 
semi-definite function on the phase space. We note that 
the Husimi distribution is not unique, but depends on 
the choice of the canonical variables (q, p) . Even for a 
specific choice of (q, p) it depends on the smearing pa- 
rameter a, as discussed below. For a two-dimensional 
quantum system, the Husimi distribution is defined as 



p H {t;qi,q 2 ,pi,p 2 ) = (z 1 ,z 2 ;a\p(t)\zi,z 2 ;a}, 



(1) 



where p{t) denotes the density operator, a is a parameter 
and the coherent state \z\, z 2 \ ct) satisfies, 
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Note that a is related to the smearing parameter A in 
@, 13 by a = h/A. The definition ([TJ ensures that the 
Husimi distribution is non-negative within all of phase 
space. Throughout this paper, the notion of p#(t;q, p) 
always implies a dependence on a, as indicated in ([1}. 
The Husimi distribution can also be obtained by Gauss 



smearing of the Wigner function. Let W be the Wigner 
function defined by: 

f°° x x 

W(t;q )P )=J d 2 x(q--|p(t)|q+-)e^- x . (4) 

The Husimi distribution is obtained by convolution of the 
Wigner distribution with a Gaussian: 



1 f°° 

p H (t; q,p) = / d 2 q'd 2 p' W(t; q', p') 

= -(q'-q) 2 A*-"(p'-p)7ft 2 
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(5) 



Since the Husimi distribution is a minimally (in the sense 
of the uncertainty principle) smeared Wigner function, it 
was proposed in Q that the Husimi distribution can be 
applied to the definition of a minimally coarse grained 
entropy, the Wehrl-Husimi entropy. In two dimensions 

/cPq cPp 
-^—^-p H (t\q,,p)lnpff(t-,ci,p). (6) 

The properties of the Wehrl-Husimi entropy are reviewed 
in [16[. Besides, Wehrl conjectured that Sjj(i) > 1 for 
any one dimensional system, where the equality holds for 
a minimum uncertainty distribution [17]. Lieb proved 
this conjecture in [l8[. We here generalize Wehrl's con- 
jecture to that of a two-dimensional system: 



S H (t) > 2, 



(7) 



where the equality holds for a minimum-uncertainty 
Husimi distribution. We confirm in Sect. HV~Al that our 
numerical results satisfy the bound ([7]). To investigate 
the time dependence of the coarse grained entropy, we 
now derive the equation of motion for the Husimi distri- 
bution. 



B. Time evolution of Husimi distribution 

In quantum mechanics, the Liouville equation 

mt) 



ih- 



dt 



(8) 



where H denotes the Hamiltonian operator, describes the 
time evolution of the density operator. One can study 
the time evolution of a quantum system by mapping the 
equation of motion of the density operator in the Hilbert 
space onto that of the corresponding density distribu- 
tion in the phase space. The Husimi equation of motion 
is obtained by subjecting both sides of eq. ([8j to the 
Husimi transform ([1]). For a one-dimensional quantum 
system, the Husimi equation of motion was first derived 
by O'Connell and Wigner [l9j . Here, we derive the the 
Husimi equation of motion for two-dimensional quantum 
system. For a single particle in two dimensions, the clas- 
sical counterpart of the Hamiltonian % reads, 
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where m is the mass of the particle and V(qi,q 2 ) is the 
potential energy. For the Hamiltonian system whose po- 
tential energy V(qi, 32) is a C°°-differentiable function of 



dt 



where A^, \ii and summed over all non-negative 

integers, with the constraints that (Ai + A2) is odd, 
(/ii — 2ki) > and (p, 2 — 2^2) > 0. When the poten- 
tial energy is of polynomial form: 

V{ qi ,q 2 ) = J2^2 a v<lUi> ( n ) 

with the coefhcients Oy and non-negative integers n\ and 
n 2 , one finds that the additional constraints (Ai + fj,±) < 
ni and (A2 + p. 2 ) _• ^2 should be applied to the sum in 

We now specialize our investigation to the Hamilto- 
nian: 

-H = ^(f 1 +pl)+\ 9 W 1 ql (12) 



It is non easy to solve the Husimi equation of mo- 
tion (fT3"|) . Before we embark on this challenge, we first 
prove the energy conservation of the Husimi function in 
Sect. Ill Cl and then solve (fT3"|) by the test-particle method 
in Sect. HD 



((71,92), we apply ((4j [5]) to flSJ, perform a series expan- 
sion of V in powers of q\ and g 2 ; and finally obtain the 
equation of motion for the Husimi distribution: 



I 

which describes a dynamical system known as Yang-Mills 
quantum mechanics [iff. The Hamiltonian in (|12[) is al- 
most globally chaotic, except for a tiny portion of the 
phase space in which stable orbits has been discovered 
[20l |2lj ]. For the potential energy in the last term of 
(IT2"1) . the order of the derivatives of V(qi, q 2 ) in (fTU)) is re- 
stricted by the relations (Ai + p,\) < 2 and (A2 + p 2 ) < 2. 
Therefore, we can rewrite the Husimi equation of motion 
([TDJ as: 



(13) 

I 

C. Energy conservation 

A coarse grained Hamiltonian, which describes energy 
conservation in the Husimi representation, was intro- 
duced by Takahashi [22Tj2^ |. who identified the quantum 
corrections to the classical Hamiltonian in powers of h 
and then constructed a conserved Hamiltonian for the 
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Husimi representation by adding these quantum correc- 
tions to the classical Hamiltonian. Explicit expressions 
for this coarse grained Hamiltonian were found for a few 
one-dimensional quantum systems [22-24]. Here wc pro- 
pose a novel derivation of the conserved coarse grained 
Hamiltonian. Our approach, which involves no approxi- 
mation, exploits the analytic properties of the transfor- 
mation between the Wigner and Husimi distributions. 

We now derive the coarse grained Hamiltonian for the 
two-dimensional Yang-Mills quantum mechanics model. 
The derivation for a one-dimensional quantum system 
is presented in Appendix |Bj Our method can be easily 
extended to the derivation of the coarse grained Hamilto- 
nian for higher-dimensional quantum systems with poly- 
nomial potentials. 

The expectation value of a Hamiltonian in the Wigner 
representation is defined as: 



/oo 
dr qiP H(q,p)W^;q,p), 
-oo 



(14) 



where T-L is the Hamiltonian, W is the Wigner function 
defined in (QJ, and 



_ rf 2 qrf 2 p 

a qp ~ (27Tft) 2 



(15) 



is the four-dimensional phase space measure. In quan- 
tum mechanics, the energy of the system is calculated as 
(H) = tv(pH). Starting from the Liouville equation flS]) 
it is straightforward to show, that d(H)/dt = 0. It is also 
easily shown [H that {%) = £[UW]. Therefore, £[UW] 
is a constant of motion under the time evolution of the 
Wigner distribution. We now apply the convolution the- 
orem to invert the transformation in ([5]) and obtain: 



£[HW] = 

where 

■Hfr(q,p) = 



dr q , p ?Mq,p)pff(*;q,p), (16) 



-L rVq'dVW,?') 
16tt 4 J_ oa 



I <i 2 ud 2 vexp 
-tu • (q' - q) - tv • (p' - p)] , (17) 



a o fi 2 n 
-u 2 H v 2 

4 4a 



and u and v are the Fourier conjugate variables to q 
and p, respectively. The expression of Hh in (IT7|) is not 
mathematically well-defined because it involves exponen- 
tially growing Gaussian functions. However, T-Lh can be 
evaluated by analytic continuation. Let £ = —a/4 and 
X] = —h 2 /(Aa). Then, we evaluate the last two integrals 
in (|17p in the analytic region where £ > and r\ > and 
obtain: 

n H (<l, p) = — -5— / d 2 q'd 2 p' W(q', p') 

1D7T ^ Tj J _ qq 

(q' - q) 2 (p' - p) 2 



x exp 



4£ 



4/] 



(18) 



Again, we evaluate the integrals in (|T5| in the analytic 
region where £ > and r\ > 0, and then we substitute 
£ = —a/4 and r\ — — fi, 2 /(4a) into its expression, thereby 
resulting in a real and finite function "H//(q, p). For ex- 
ample, by substituting (H~2l) into (H"8|) and evaluating (fT8|) 
according to the above procedure, we obtain: 



H*(q,p) = ^ 
1 



(Pi 



2\ i 1 2 2 2 
P2j + 2^ ^1^2 



-<?-a (g 2 + g 2 ) 

1 2 2 ^ 2 
o3 a 
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The analytic function %ij(q, p) in (|T^|) is the coarse 
grained Hamiltonian for the Yang-Mills quantum system 
whose conventional Hamiltonian is defined in (IT2|) . We 
now define the expectation value of the energy in the 
Husimi representation as: 



/oo 
dT^p Ktf(q,p)pff(i;q,p), 
-OO 



(20) 



where Hh (q> p) is the coarse grained Hamiltonian de- 
fined in flTHJl. Using eqs. ([T^ 1 H5 1 [2D | . it is straightforward 
to prove by explicit calculation that 



d£[n H p H ] 
dt 



= 0. 



(21) 



Thus, £[HhPh] is a constant of motion for the Husimi 
equation of motion (|13[) and can be identified as the total 
energy of the system. In Sect. lIV"^! we verify numerically 
that £ [HhPh] is a constant of motion. 



III. TEST PARTICLE METHOD 

The numerical solution of the Husimi equation of mo- 
tion for one-dimensional quantum systems has been in- 
vestigated, e.g., in [26], |27J. Because our goal is to ap- 
ply the Husimi representation to quantum systems in 
two or more dimensions, we need a method that is ca- 
pable of providing solutions to the Husimi equation of 
motion for higher-dimensional systems. As a practical 
approach to this problem, we here adopt the test-particle 
method. This method was previously applied by Heller 
28], who assumed that the wave function is a superpo- 
sition of frozen Gaussian wave packets. The test-particle 
method was also used to describe the time evolution of 
the Husimi function of one-dimensional quantum systems 
by Lopez, Martens and Donoso [27[. Manipulating the 
Husimi equation of motion algebraically, these authors 
obtained the equations of motion for the test particles. 
The equations of motion for test particles obtained in this 
manner exhibit a nonlinear dependence on the Husimi 
distribution. However, we note that the true equation of 
motion for the Husimi distribution is a linear partial dif- 
ferential equation, which encodes the superposition prin- 
ciple for quantum states. The nonlinear dependence of 
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the equations of motion for the test-particles represent- 
ing the Husimi distribution in [27| implies a violation of 
this principle. We note that the superposition princi- 
ple is crucial to our investigation. To study the entropy 
production of the Yang-Mills quantum system and the 
approach to thermal equilibrium, we need to consider 
highly excited states of the system, whose energies form 
a quasi-continuum. Thus, the time evolution of the sys- 
tem is described by the superposition of eigenstates with 
almost the same energy. When the superposition prin- 
ciple is violated, we cannot expect to describe the time 
evolution of such states correctly. 

Therefore, we here apply the test-particle method in 
a way that respects the superposition principle. Instead 
of adopting the strategy proposed in [27| , we obtain the 
equations of motion for the test particles by taking the 
first few moments on the Husimi equation of motion. 
This approach preserves the superposition principle for 
solutions of the Husimi equation of motion. In Sect. lIII Al 
we derive the equations of motion for the test particles, 
obtain the uncertainty relation for Husimi distribution, 
and prove that the energy conservation holds for each 
individual test particle. In Sect. IIIIB1 we describe the 
method by which we choose the initial conditions for the 



Husimi equation of motion. In Sect. IIII C[ we discuss 
additional approximations that we use for the Gaussian 
test functions. 



A. Equations of motion for the test particles 

Now we briefly describe the test-particle method. Our 
goal is to solve the Husimi equation of motion in (| 13[) and 
obtain the time dependence of the Husimi distribution. 
As stated before, the Husimi distribution is a density 
distribution on the phase space, and it is positive semi- 
definite for all times. Therefore, we can approximate the 
time-dependent Husimi distribution by the superposition 
of a sufficiently large number N of Gaussian functions, 
whose centers can be considered as the (time-dependent) 
positions and momenta of N "test particles" . 

For these Gaussian functions, we assume that we can 
neglect all correlations between qi and q 2l between p\ 
and p 2 , between qi and P2, and between q 2 and p\. Un- 
der these assumptions, the Husimi distribution can be 
written as 



Pff(*;q,p) = — ^V^W ex P 



i=l 



1 



x exp 



-A lpl (t)( Pl -p\(t)) 



4 ft (*) fa -£(*)) 
1 



2 ?2<?2 



(*) fa-titty 



2 C P2P2 



(t){ P 2-P l 2 (t)Y 



x exp [-c qxpi (t) ( qi - ql(t)) ( Pl - p\(t)) - 4 P2 (t) (q 2 - 4(t)) (pa - p 2 (t))} 



(22) 



In order to satisfy the normalization condition for the 
Husimi distribution: 



/oo 
dr q , p p H (q,p;t) = 1, 
-OO 



we normalize each Gaussian according to: 

N*(t) = A*(t)A*(i), 
where we introduced the abbreviations: 



A'i(i) 



49i(*)4px(*)-(4px(*)) 



m) = 4»(*)4»(*)-(4p a (*))' 



(23) 

(24) 

(25) 
(26) 



We require that N l (t) > for all times. The assumption 
of setting c\ iqa (t) = 4 P2 (i) = c 9iP2 W = c 92 P i W = in 
(f22j) is motivated by the fact that c l qipi (t) and c q2P2 (t) en- 
code the dominant correlations induced by the dynamics. 
For purposes further down, we have examined numeri- 
cally that even when setting c* = 2 (t) = for 



all times, the correlations between qi and p\ and between 
q 2 and p 2 are produced by the ensemble of Gaussians 
as time evolves, by virtue of the contribution of a large 
number of test functions. Therefore, the ansatz in (|22|) 
is justified. 

Owing to (|22|) . the solution to the Husimi equation of 
motion will depend on the chosen particle number TV, and 
so will the Wehrl-Husimi entropy. In the limit N — > 00 
we expect both, the Husimi distribution and the Wehrl- 
Husimi entropy, to approach values that are independent 
of the test particle approximation scheme. We will con- 
firm this expectation in Sect. II V CI by investigating the 
particle number dependence of our numerical result for 
the Wehrl-Husimi entropy. 

The main task for us is to determine the optimal solu- 
tions for the time-dependent variables q\(t), q\(f), P\(t), 



C 9l9lW' 4292 (0' C PlPlW> 



,(*) 



' 9lPl 



(t), and 



c l q2P2 (t). In other words, instead of directly solving (fl3|) . 
we seek a system of the equations of motion for the 
ten time-dependent variables. This goal can be achieved 
by evaluating the moments on both sides of the Husimi 
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equation of motion. The resulting equations constitute 
a system of ordinary differential equations for the ten 
time-dependent variables of each test particle labeled by 
i = 1, 2, N. Overall, we thus have to solve 10 N equa- 
tions of motion. These can be grouped into N indepen- 
dent systems of ten coupled differential equations, each 
of which can be solved separately. 

Generally, the moment of a function f(t; q, p) with re- 
spect to a weight function iw(q, p) is defined as, 



dr q , p [w(q,p)/(i;q,p)] 



(27) 



Therefore, after we apply the ten moments I qi , I q2 , I Pl , 



T I 2 



I qipi and Iq 2P2 to the Husimi equa- 



tion of motion (|13l) , we obtain ten equations of motions 
for each test particle i for the ten variables representing 
the location in phase space and width of each test par- 
ticle. In eqs. (|28M31I) . we present the equations obtained 
from the first moments I qi , I q2 , I Pl and I P2 of (| 13[) associ- 
ated with the location of the test particle. The equations 
for the evolution of the test particle widths, obtained 
from the second moments I q 2 , I p 2 , I p 2 , I qiPl and I q2P2 of 
([TB")) are presented in eqs. (|A1IIA6[) of the Appendix |A"1 
The equations for the first moments of (fTB"j) are: 



+ 


dV 






dqi 


q*(t) 


\ 






V a| (t) 




dV 




+ 




q*(t) 




c pipi W 




V Aj (t) 



d 3 V 



dqidq. 



= 0, (28) 
= 0, (29) 

= 0, (30) 



q'(t) 



d 3 V 



0, (31) 



q<(t) 



dq\dq 2 



where Ai(t) and A l 2 (t) are defined in <(25]) and ((26j) , re- 
spectively. The subscript q J (i) in the partial derivatives 
of the potential energy V(q%,q2) in (f30"ll3Tj) denotes that 
the partial derivatives are evaluated at (91,92) — q l (^)i 
where 



(32) 



(«) = $(*),§£(*)) 



Instead of solving the Husimi equation of motion f|13[) . 
we now solve (|28H31[) and (|A1IIA6P for each test particle 
i = l,2,...,iV and then construct the Husimi distribu- 
tion by superposition. These test particle equations of 
motion can be solved numerically by applying the Runge- 
Kutta method when proper initial conditions are given. 
The method of choosing the initial conditions will be dis- 
cussed in Sect. IIIIB1 

To ensure the existence of the solutions, we need to 
confirm that eqs. (|AHIA6p are nonsingular. We write 



the system of differential equations (IA1IIA6I) in the form 
Av = b, where v and b are column vectors and 



The system of equations would be singular if det A 
which implies, 

A*(i)A*(t)=0. 



(33) 
= 0, 

(34) 



This condition is equivalent to N l (t) = 0. Equation (|3"4")l 
violates the constraint that N l (t) > 0; therefore, (|28H31[) 
and (|AHA6I) are never singular. 

The uncertainty relation for the Husimi distribution 
for one-dimensional quantum systems has been derived 
in, e. <?., [1^1 . Here we generalize their result to the case of 
two dimensions. The uncertainty relation for the Husimi 
distribution pu (t; qi,q2,Pi,P2) reads: 



(A<&)h (Apj) H > h, 



(35) 



where 



(Aft)* 

for j = 1,2 with 

{<1j)h 



xp H (t;q,p)} 



(Pi 



AcK*;q,p)] , 



/oo 
-OO 

/OO 
-OO 



(36) 
(37) 

(38) 
(39) 



We emphasize that the uncertainty relation ([33)) does 
not serve as an additional constraint when we solve the 
Husimi equation of motion (|13[) . As long as the initial 
condition pa (0; q\, 92,^1,^2) satisfies (|35[) . the solution 
to the Husimi equation of motion satisfies the uncertainty 
relation (1351) for all times. This results from the fact that 
the quantum effect is encoded in the Husimi equation of 
motion itself. 



B. Initial conditions 

In order to solve the equations of motions (|28H31l 
IAHIA6I) . we need to assign initial conditions for the 
Husimi distribution at t — 0. We next describe 
the method we use to assign the initial conditions, 
{<Zi(0), q\{fyiP\ (0),P 2 (0)} an d the initial widths for each 
test particle i. Our goal is to assign initial conditions 
so that the initial Husimi distribution satisfies the four 
conditions at t — 0: (i) /Ojj(0;q, p) > 0, («) the normal- 
ization condition in (|23[) . (Hi) the uncertainty relation in 
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(l35f . and (iv) the relation between moments: 

/OO />00 
dr q , pjO tf(0;q,p) > / dr q , p [p ff (0;q,p)] 2 .(40) 
-oo J —oo 

Our strategy is as follows. First of all, we formally write 
([22]) as: 

1 N 

M*;q,p) = ^^^(q-q l W,P-p 4 W), (41) 

i=l 

where if denotes the Gaussian function for each test par- 
ticle. For t = 0, the Husimi distribution (|4"Tj) can be 

expressed as 

/oo 
dTq'.p' if (q - q', p - p') 
- OO 

X0(q',p'), (42) 

where (j> denotes the distribution of the test particle lo- 
cations in the phase space. We abbreviate the phase 
space variables for clarity: x = (Qi> Q2,Pi,P2) and x' = 
(si) <?2> Pi iP^)- Owing to the four conditions (i)-(iv) 
stated above, we choose the Husimi distribution at t = 
to be a Gaussian distribution: 



ph(o ;x ) = n 2 




x exp 



a=l 



(43) 



where jfj and p a K for a = 1, . . . , 4 are to be determined. 
In (|43j) we do not assume any correlation between posi- 
tion and momentum locations at t = 0, implying that we 
initially set c4 ipi (0) = c l X2P2 (0) = for i = 1, . . . , N in 



The main idea of choosing initial conditions is that, 
according to (|42p . we can represent the initial Husimi 
distribution (1431) to be the sum of Gaussian test functions 
by randomly assigning {qi(0), ^(0),pl(0),p|(0)} for i 
1 . N according to the distribution (f>. Our remaining 
tasks are then to determine the parameters in (|43l) and 
to obtain the functional forms for K and </>. In (|4"3"]) . p a H 
can be assigned freely by choice, but the 7^ are subject 
to the conditions (Hi) and (iv). Substituting (1431) into the 
conditions (Hi) and (iv), expressed by eqs. (|35|) and (J40J) , 
respectively, we obtain from (Hi): 



ri(-y&) 



-1/2 



> h 2 



(44) 



and from (iv): 



n (7&r 1/a > ^ 2 /4. 



(45) 



Since eq. (l44l) is the stronger of the two conditions, we 
adopt it as the constraint for the initial Husimi distribu- 
tion. To represent ph(0, x) m (|43| . we chose the follow- 
ing functional forms for K and tf> at t = 0: 

, 4 xV2 

^(x-x') = ^ 



x exp 



1 4 

-^E^(x a -x' a r 



,(46) 



and 



1/2 



hx) = a 2 



x exp 



a=l 



(47) 



This choice implies that we represent the initial Husimi 
distribution as the convolution of test particle Gaussian 
functions K and a Gaussian distribution <f> of test particle 
locations in phase space. In (|4T|) at t = 0, pn is denoted 
as the sum of Gaussian functions, each of which may 
possess distinct widths. However, when we choose to 
express (|4Tj) at t — in terms of the convolution of K 
and </>, we no longer have the flexibility to assign different 
widths for each individual Gaussian. Instead, for K in 
(|4*21 HSJ) we should assign: 



Ik = c 9igi (0), lK = c q2q2 (0), 



IK 



1 in 



Pl (0), 



IK 



U P2P2 



(0), 



(48) 



where the suppression of the label i implies that all test 
particles possess the same width at t = 0. 

It is advantageous to use the convolution of K and </> in 
(|4*21 to represent ph because the parameters in (|4"51 |4U1 
|4T|) can be related to satisfy the constraint imposed by 
the uncertainty condition, as described below. In (|47p . 
/Lt^ denotes the location of the center of the distribution 
of loci of the test particles in the phase space. According 
to (g21 |H |H ST]), it is clear that the center of the dis- 
tribution of loci of test particles must coincide with the 
center of the initial Husimi distribution. We thus must 



assign 



(49) 



where [i a H are selected by choice. Moreover, since the 7^ 
are subject to the constraint (|44p . we obtain relations 
between 7^, j K and 7^, which allow us to determine j K 
and 7^. By applying the convolution theorem to (|42j) . we 
obtain the following relations: 



1 

1% 



1 

7& 



1 



(50) 



for a = 1, . . . , 4. Once we select the values of 7^ based 
on (|4~4"|) , we must determine 7^- and 7? according to (j5T 
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Furthermore, owing to (|48|) . the choice of 7^ is subject 
to the constraints 

1k>1h for a =1,..., 4. (51) 

Furthermore, 7^ must be assigned in the domain where 
the solutions of (j28I31|) and (|A1IIA6I) are stable. We dis- 
cuss our choice of initial conditions in more detail in 
Sect. UVA"! 

The number N of test particles plays a crucial role for 
the accuracy of numerical results. If we set N = 1 in 
(|4Tj) . we find that pn = K, and thus — Ik- This spe- 
cial case is called the single- particle ansatz. In general, 
the single particle ansatz is insufficient as representation 
of Ph (t; qi, <Z2jPijJ?2)j because the Husimi distribution 
will not retain a Gaussian shape for all times, even if we 
initialize it as a Gaussian at t = 0. 

As a specific example, we present and compare the 
solutions of the Husimi equation of motion in one di- 
mension in Fig. Q3] in Appendix [Bj Figure [H] shows the 
difference between the solution pH(t;q,p) for the single 
particle ansatz [panels (a) and (b)] and for the many- 
particle ansatz [panels (c) and (d)], for the same Hamil- 
tonian defined in eqs. (|B2[) . The initial conditions are 
also discussed in Appendix [Bl From Fig.fTJl it is obvious 
that the single-particle ansatz is insufficient in represent- 
ing the solution pn(t;q,p) for t > 0. We conclude that 
we need a sufficiently large test-particles number N in 
(|4"Tj) to represent the evolution of the Husimi distribu- 
tion. We discuss the test-particle number dependence of 
our numerical results in Sect. IIVI 

C. Fixed-width ansatz 

Once the initial conditions are obtained, the numeri- 
cal solutions to eqs. (|28M31I IA1MA6|) can be obtained by 
the Runge-Kutta method. These equations can be dra- 
matically simplified by fixing the Gaussian widths in our 
ansatz ([22)1 for the Husimi distribution. The precise defi- 
nition of the fixed-width ansatz reads as follows: For each 
particle i, 





= c 9l<3l(0)> 




— C <J2<?2 




= c pipi(0), 


<£«»(*) 




4px(*) 


— c qipi(0), 




= C <?2P2(0)j 



where c gigi (0), c q2q2 (0), c plPl (0), c P2P2 (0), c 9lPl (0), and 
c ?2P2 (0) are chosen to be the same for all i. 

In the variable- width ansatz, we solve the ordinary 
differential equations (|28H31I IAHIA6[) simultaneously for 
each test particle i. In the fixed-width ansatz, we fix 
the values of c qiqi (t), c q2q2 (t), c* m (t), c p2P2 (t), c qipi (t), 
and c qipi (t) to be constant for t > 0. Therefore, in the 
fixed-width ansatz, eqs. (|AHIA6I) cannot be satisfied, and 
eqs. (|28H3ip are the only equations of motion for each 
test particle i. We apply the fixed- width ansatz because 
(|28H3ip are obtained from the first moments of (fT"3"|) and 



thus serve as the leading contribution to (|T3|) . From a 
physical viewpoint, equations (|28H3ip determine the "lo- 
cations" of test particles in the phase space as functions of 
time, while eqs. (|AHIA6|) govern the time- varying widths 
of each test-particle Gaussian. In Sect. IIVI we evaluate all 
of the numerical results based on the fixed-with ansatz 
in (52). 

The conservation of energy is not only true for pn , as 
shown in Sect. Ill CI but also holds for each individual 
test particle. We now prove the conservation of energy 
for each individual test particle in the fixed-width ansatz. 
The proof can be easily generalized to the case of variable 
widths. In the fixed-width ansatz, the test-particle space 
is spanned by the test-particle positions and momenta 
(q, p). We define a function Hh in the test-particle space 
as follows: 

/oo 
dTq.p H H (q, p) K (q - q, p - p) , 
-OO 

(53) 

where Hh denotes the coarse-grained Hamiltonian de- 
fined in Sect. Ill CI and K is defined in (|4"Tj) . We note 
that the functional form of K is independent of the test- 
particle label i. With the help of (J2TJ) and ([55]). it is 
straightforward to show that 

8u H mwm = ^ (54) 

where i = 1, N. In view of fUJ, U H (cf(*), P*(*)) can 
be identified as the energy of an individual test particle 
i. Due to (|54|) . the histogram of test particle energies 
Hh (q z (^)j P z (^)) remains unaltered at all times. We ap- 
ply this result to the numerical calculation in Sect. IIVI 

Before we end this Section, a general consideration 
is in order. In principle, any smooth, positive definite, 
normalizable function on the phase space can be repre- 
sented to any desired precision by a sufficient number of 
sufficiently narrow Gaussian functions with fixed width. 
However, it is important to keep in mind that these con- 
ditions are not satisfied, in general, by the Wigner func- 
tion or the classical phase space distribution of a chaotic 
dynamical system. The Wigner function is in general 
not positive definite, and the classical phase space dis- 
tribution does not remain smooth for an arbitrary initial 
condition. The presence of exponentially contracting di- 
rections in phase space ensure that, over time, the clas- 
sical phase space distribution will develop structure on 
exponentially small scales, which cannot be described by 
superposition of fixed-width Gaussian functions. 

The Husimi transform of the Wigner function cures 
both problems. It removes regions of negative values 
from the quantum phase space distribution, and its re- 
spect for the uncertainty relation ensures that the phase 
space distribution remains smooth on the scale set by h 
and the smearing parameter a. As a result, the fixed- 
width Gaussian ansatz will always be able to represent 
the Husimi distribution and track its evolution faithfully 
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over time, if a sufficiently large number of sufficiently nar- 
row Gaussian test functions is employed. On one hand, 
the width of Gaussian test functions cannot be larger 
than the width of the initial Husimi distribution so that 
the Gaussian test functions can represent pn faithfully, 
as indicated in (f5fjl . On the other hand, the width of 
Gaussian test functions must not be too narrow in order 
to ensure that the solutions of (|28M31 [) are stable. We do 
not attempt to give a rigorous proof of these assertion 
here, but content ourselves with the heuristic argument 
presented above. We will explore the convergence of or 
numerical solution for the fixed-width ansatz for large 
values of N at the end of the next Section. 



IV. NUMERICAL RESULTS 

We now present our numerical results. Throughout 
our calculations, we have used the fixed- width ansatz as 
described in Sect. IIII CI In Sect. IIV Al we present the 
numerical results for the evolution of the Husimi dis- 
tribution and the Wehrl-Husimi entropy of the Yang- 
Mills quantum system using N — 1000 test particles. 
In Sect. IIV B| we obtain the Lyapunov exponents, the 
average Kolmogorov-Sina'i entropy and the logarithmic 
breaking time for Yang-Mills quantum mechanics. In 
Sect. IIV C| we compare the Wehrl-Husimi entropies for 
N = 1000 and N = 3000 test particles and explore the 
test particle number dependence of the saturation value 
of the Wehrl-Husimi entropy. In Sect. IIV D) we obtain the 
partition function and entropy for the canonical ensem- 
ble. Then, in Sect. IIV El we evaluate the microcanonical 
distribution and entropy, and we compare the saturated 
Wehrl-Husimi entropy to the microcanonical and canon- 
ical entropies. 



A. Husimi distribution and Wehrl-Husimi entropy 

For our numerical calculations, we fix the parameters 
m = g = a = h= Xixi (|I3p . Initially, we set the number 
of test particles to N — 1000. We choose a minimum 
uncertainty initial Husimi distribution (|43[) by setting: 



for a = 1 , 



(55) 



which satisfies the constraint (|44|) . Besides, in (|43p we 
choose 



2 



0, p% =fi* H = 10. 
Owing to (l4T)l 151))) , we then have 

v\ = tfj> = 0, [i\ = [i\ = 10. 



(56) 



(57) 



For a fixed-width ansatz, the solutions of (|28ll31j) are sta- 
ble under the following constraint: 



l(0) 



,(0) 



(0)c g292 (0) 



>a, 



(58) 




8 10 



FIG. 1: Conservation of the coarse grained energy (|2U|I during 
time evolution of the Husimi distribution. This shows that a 
state with energy £ [HhPh] ~ 100.707 for t = remains at 
the same energy for t > 0, with relative precision better than 
10" 4 up to t = 10. p H is obtained from @IJ with N = 1000 
fixed-width test particles. 



which can be confirmed by a linear stability analysis. 
Besides, we set c 9lPl (0) = c 92P2 (0) = according to 
Sect. iHIBl Thus, due to (p) 155)1. our choices of 7^ 
and 7^- are constrained by: 



IkIk 



> a. 



(59) 



In summary, our choice of 7^- is restricted by the two 
constraints ((511 159")) together with the settings ((55)) and 
a = 1. In view of the discussion in Sect. IIII Bt we satisfy 
these constraints by the choice 



IK 



7^ = 3, (a=l,...,4). 



(60) 



As described in Sect. IIII Bt we randomly generate test 
particle locations {q[ (0), ^(0),p{ (0) ,p 2 (0)} f° r * = 
1 , . . . , N according to <fi in ((47)) , with parameters given 
by ((571 IBT))) . For the fixed- width ansatz with the initial 
conditions (|60)) . we solve (j28H31l) for each test particle i 
and repeat the procedure for i = 1,2, N. 

Using cqs. (fl"9l 120)) where pn is obtained from ((4T)) 
with TV = 1000 fixed-width test particles, we verify nu- 
merically that £ \HhPh] is a constant of motion. This is 
illustrated in Fig.[l] which shows that a state with initital 
energy £ \HhPh] = 100.707 remains at the same energy 
with relative precision better than 10~ 4 up to t = 10. 
Since the initial " locations" of test particles in the phase 
space are generated randomly according to <f> in (|47[) , dif- 
ferent sets of {q l (0),p l (0)} generated by different runs 
of the computer program may result in differences of 
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FIG. 2: Energy histogram for N = 1000 test particles at t = 0. 
e denotes the test-particle energy, which is defined in (|6ip . and 
the labels on the vertical axis denote test-particle numbers. 
A normal distribution nTp(e) is used to fit the histogram. A, 
fi and a are the fit parameters for riTp(e), which are defined 
in (I64[) . The values for the fit parameters are shown in the 
plot. 



£ [HhPh] at t = of less than 0.5 percent. Thus, for 
any set of initial locations for N — 1000 test particles, 
the energy of the state at t = is £ [HhPh] = 100.6±0.5. 

The energies of individual test particles can be studied 
by the following method. We denote the test-particle 
energy variable e as 



£ = Hh (q,p) 



(61) 



where Hh (q, p) is defined in (|53)) . Because we choose 
the fixed-width Gaussian K with the parameters 7^- in 
(|60l) and set m = g = a = K=l, we obtain 



Wh (<Zl,«2,Pl,P2) 



(Pi 



1 

2 

1 /-a 
12^ 



91) 



-2 -2 



13 

72 



^- (62) 



The energy for an individual test particle is denoted as i 
e* = Kff (q* Owing to gl]), the energy of the 
state is the average energy of the test particles: 



£ [HhPh] 



1 N 
N ^ 

i=l 



provided that N is sufficiently large, 
the energy histogram at t = for N — 
which wc fit to a normal distribution: 



(63) 



In Fig. [2] we plot 
1000 test particles, 



^tp (e) = ^4exp 



1 

2^2 



(64) 



The values of the fit parameters A, /i and a are listed 
in Fig. [2] for N = 1000. We note that the histogram of 
test particle energies remains unaltered as time evolves, 
as shown in Sect. MI CI 

To visualize the Husimi distribution as a function of 
time, it is useful to project the distribution either onto 
the two-dimensional position space (91,92) or onto mo- 
mentum space {p\,P2) by integrating out the remaining 
two variables. To this end, we define the following two 
distribution functions: 



F q (t;q 1 ,q 2 ) 



/OO 
-00 



2nh 2 
N 



AiA 2 



c PlPl c P2P2 



exp 



Ai 

2cp 1 



A, 



2c 



P2P2 



(65) 



/OO 
dqidq 2 p H (t;q 1 ,q 2 ,pi,p 2 ) 
-CO 



-00 

2itK- 



E 



AiA 2 



■ exp 



Ai 



2c 



qiqi 



(Pi-P\(t)f 



2c, 



[P2-m)Y 



qiq-i 



(66) 



We can conveniently visualize the evolution of the 
Husimi distribution pn(t; qi, q 2 ,p\,p 2 ) by showing con- 
tour plots of the two-dimensional projections F q (t; q\ , q 2 ) 



and F p (t;pi,p 2 ). Figure [3] shows F q and F p side by side 
at times t = 0, t = 2, and t — 10, respectively. At 
the initial time, F q (0; qi, q 2 ) is chosen as a Gaussian dis- 
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tribution centered around the origin in position space, 
while F p (0;pi,p2) is a Gaussian function centered around 
(p 1 ,p 2 ) = (10, 10). The projected initial distributions are 
shown panels (a) and (d) of Fig. [31 As shown next in pan- 
els (b) and (e) of Fig. [3j F q and F p at t — 2 are beginning 
split into distinct clusters. This behavior is caused by the 
fact that test particles bounce off the equipotential curves 
defined by e = FLh (q, 0). 

Closer inspection of the time evolution of F g (t;qi,q 2 ) 
and F p (t;pi,p2) reveals that gross features of the Husimi 
distribution pnit', q, p) remain approximately unchanged 
for t > 6. To wit, the panels (c) and (f) of Fig. [3J 
presenting F q and F p at t = 10, show that the con- 
tours of F q (10; qi, q 2 ) follow equipotential lines, while the 
contours of F p (10;pi,p2) are shaped as concentric cir- 
cles, i. e., lines of constant kinetic energy. The time 
evolution of F q demonstrates that test particles start- 
ing from their initial positions localized around the ori- 
gin in position space ((71,(72) eventually spread all over 
the region enclosed by the equipotential curves defined 
by e = Hh (q, 0). This behavior is a result of the fact 
that the Yang-Mills quantum system is chaotic implying 
a strong sensitivity of test particle trajectories on their 
initial conditions. 

The Wehrl-Husimi entropy Sn(t) defined in © is the 
coarse grained entropy of the quantum system. The nu- 
merical evaluation of the four-dimensional integral in the 
definition ([6]) of the entropy Sn{t) is nontrivial because 
the upper (lower) limits of the integral in each dimension 
are infinite and the width of each test particle Gaussian 
is tiny. Therefore, we use the following method to eval- 
uate the integrals efficiently. For each discretized time 
step tk, we find the largest absolute values of the test 
particle positions and momenta. Since each Gaussian is 
narrow and the Husimi distribution is nearly zero out- 
side the regions of support of the test particles, we can 
assign ±(max^ |<?j(ifc)| + b) as the limits of integration 
over the variable q%. We choose b = 6(7] f )~ 1 / 2 to ensure 
that the Gaussians of all test particles are fully covered 
by the integration range within our numerical accuracy. 
Similar limits are assigned to the integrations over q 2 , p\ , 
and P2 , respectively. These integration limits ensure that 
the integrals run over the whole domain of phase space 
where the Husimi distribution has support. We verify 
the accuracy of Simpson's rule by evaluating the normal- 
ization for pnit', q, p) for various time t. We find that the 
numerical results coincide with (I23[) within errors of less 
than 0.3%. We then perform the numerical integration 
by Simpson's rule. 

Our results for the Wehrl-Husimi entropy Sjj(t) for 
N = 1000 test particles are shown in Fig. |4j We evalu- 
ate Sji(i) for Yang-Mills quantum mechanics (YMQM) 
and for the harmonic oscillator (HO), for comparison. 
The Hamiltonian for YMQM is given in (|12p. while the 
Hamiltonian for HO is: 

n = ^-(pl+pl) + 1 -v 2 ( q j + ql), (67) 
2m 2 

where we set m = v = 1. We remind the reader that 



initially /9jj(0) is chosen as a minimum uncertainty dis- 
tribution satisfying the constraints fHl [55)1 with the to- 
tal number of test particles N — 1000. We assign the 
same initial condition both for YMQM and HO, and we 
compare the difference in their Wehrl-Husimi entropies 
as time evolves. Figure |4] shows that Sh{0) » 2.0, and 
Sh(0) > 2 for t > for YMQM, in agreement with the 
conjecture For late times, Fig. @] reveals that Sjj(t) 
for YMQM saturates to 7.7 for t > 6.5. In order to find 
the characteristic time for the growth of the entropy, we 
fit Sji(i) for YMQM to the parametric form: 

5st(t) = «o-«iexp(-t/r), (68) 

where sq, s± and r are fit parameters. The fit shown as a 
dash-dotted line in Fig. 0] corresponds to the parameters 
so ~ 7.7, si « 6.0 and r w 1.9. On the other hand, Sn(t) 
for HO starts from Sh(0) ~ 2.0 and then remains at 2.0 
for all times. 

In Fig. HI we note that the coarse grained entropy 
does not increase continuously as time evolves. This fact 
can be interpreted in the framework of Zwanzig's formal- 
ism for the time evolution of "relevant" density operator 
[29l |30T |. In Zwanzig's formalism, one defines the rele- 
vant density operator as pn(t) = Pp(t), where P denotes 
the projection operator. The transition of the density 
operator p(t) — > Pn(t) and of corresponding entropies 
S[p(t)] — > S[pa(t)] is referred to as generalized coarse 

graining [3(J[3l|. a Pply m g P to ©, one obtains the 
equation for time evolution of pn(t). The non-Mar kovian 
part of this equation reads: 

^M = -J t dt'G(t')p R (t-t% (69) 

where G denotes the so-called memory kernel [29j - [3ll | . 
It can be shown that dS[pfj(t)]/dt receives contributions 
from the non-Markovian term indicated in (l69l) . There- 
fore, S[pn(t)] in general does not increase monotonically 
as a function of time. The Husimi equation of motion 
in (|10[) contains a similar memory effect. Therefore, in 
Fig.|4]the coarse grained entropy Sn(t) does not increase 
continuously as time evolves, and the second law of ther- 
modynamics holds only in a time averaged sense [30| . 

B. Lyapunov exponents 

Since the classical system corresponding of YMQM 
is almost chaotic, we evaluate the average Kolmogorov- 
Sina'i (KS) entropy for this system. For a two dimensional 
system, the KS entropy is defined as: 

4 

fcKS = 5>i*( A i)> (70) 
3=1 

where Aj's are the Lyapunov exponents (LE). To obtain 
the full spectrum of the LEs, we utilize the following pro- 
cedure. First, we divide a large time interval, from t = 




FIG. 3: Two-dimensional projections of the Husimi distribution on position space F q (t; qi, 92) at times (a) t = 0, (b) t = 2 and 
(c) t — 10, and on momentum space -F p (i;pi,p 2 ) at times (d) t = 0, (e) £ = 2 and (f) t — 10. The number of test particles is 
N = 1000. 
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FIG. 4: The time evolution of the Wehrl-Husimi entropy 
Suit) for Yang-Mills quantum mechanics (YMQM), the fit 
function Sfit(t) for the Wehrl-Husimi entropy, and Sn{t) for 
the harmonic oscillator (HO). We set the same initial con- 
dition at t = both for YMQM and HO. The figure shows 
that S H (t) for YMQM starts from Sh(0) ~ 2.0 and saturates 
to 7.6 for t > 6.5, while Sn(t) for HO remains at 2.0 for all 
times. The fit parameters for Sat(t) are listed in the figure. 



to t = t max , into a number of slices. Each time slice 
is labeled by its final time tk, where k = 1, 2, fc max . 
Let y?{t) = {^(t),^(t),Pi(t),p2(t)) denote the position 
of test particle i in phase space. At t — 0, we per- 
form four orthogonal perturbations on the initial con- 
dition: 7r](0) = x'(0) + for j = 1,...,4, where e/s 
are orthonormal vectors and we set e = 10 -4 . For each 
time slice t € t^l, we solve eqs. (128H31I) and obtain 



i,£fc], we solve eqs. (|28H31|) 
one reference trajectory % z (i) and four modified trajec- 
tories TTj(t), where j = 1,...,4. Define the four devia- 
tion vectors: 8j(t) = 7r*(i) — X 4 (£). After obtaining the 
four deviations <5j(ifc), we orthogonalize these four vec- 
tors and rescale their lengths back to e. We store the 
four rescaling factors r*(tk) for each j and k, and we 
repeat the above procedures for the representative test 
particles i — l,...,N Tep , where A rop < N. For the case 
of N = 1000, we choose iV re p = 100. Besides, we set 



tk = 2fc and t D 



100, and therefore k B 



nally, we obtain the full Lyapunov spectrum: 



l l 

rep - . f'ma 



In 



n 

fc=i 



50. Fi- 



(71) 



where j = 1, 
YMQM are: 

Ai = 
A 3 = 



, 4. The numerical values of the LEs for 



1.216, A 2 = 2.344 x 10 
-2.349 x 10~ 2 , A 4 



-2 

1.223. 



(72) 



If we take the classical limit h — > and a — > for and 
repeat the above procedure, we obtain the LEs for the 
regular classical equations of motion without the quan- 
tum (Husimi) corrections: 



Ai 

A, 



1.283, A 
-1.629 x 10 



2 _ 

-2 



1.599 x 10 , 
X% = -1.287. 



(73) 



From (|72j) and (f73|) . we observe that classical solutions 
conserve the energy and the phase space better than the 
quantum solutions. By (|70l l72j) . we obtain the average 
KS entropy for YMQM: h KS w 1.24. 

In addition, we calculate the logarithmic breaking time 
for YMQM, which is defined as [321131}: 



Hi 



(74) 



where / is the characteristic action and A is the charac- 
teristic Lyapunov exponent. We set A = /iks for YMQM. 
We utilize two methods for obtaining the action /. One of 
these is to obtain / from the classical dynamical variables 

(q,p): 



I = 



(75) 



The integration is taken over the curve C constrained by 
H = E, where % is defined in (fl~2"|) and E denotes the 
classical energy of the system. If we consider the case 
where a classical particle moves along the line q± = q2 in 
the position space and is subject to the potential energy 
|<?i<72j we obtain the period of motion of this classical 
particle: 



T 



dq 



E 



h 4 



(76) 



where q = q\ = qi and (?max = (2E) 1 / 4 . In the following 
numerical calculation, we set E = 100. Considering the 
periodic motion of this particle, we obtain by (I74H76[) 
that I = 263, T = 1.97 and m w 4.5. Alternatively, we 
evaluate the action by the integrating along test-particle 
trajectories obtained by (|28H31[) , which are the Husimi 
(quantum) equations of motion in the fixed- width ansatz. 
Thus the action is: 



I = 



N x 

1W dtp*®.?®, (77) 

,- i */ 



where T defined in (|76|) . In (f77|) . we estimate the time 
interval by the period of a classical particle moving along 
Qi = 92 in the position space and is subject to the po- 
tential energy \q\q\. By ([711 I77|l, we obtain I = 267 
and Th f=a 4.5 in excellent agreement with the result of 
the first method. Moreover, comparing r% to r defined 
in (I68p . we conclude that t% and r are in the same order 
of magnitude, and > r. 
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choices of N, we obtain 
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(78) 
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FIG. 5: Energy histograms of the test particles at t = 0. 
The total numbers of test particles are iV = 1000, N = 3000, 
and N = 8000. e denotes the test-particle energy, which is 
defined in (|61[) , and the labels on the vertical axis denote test- 
particle numbers. The initial locations of the test particles in 
the phase space are generated according to the normal dis- 
tribution <f) defined in (|47[) with the parameters given in (|57l 
I60|l . In this plot, we show that fi and a are independent of 
N, notwithstanding small fluctuations. By fitting the energy 
histograms for various choices of N, we obtain n = 100.6 and 
a — 8, with fluctuations less than 0.5% and 5% respectively. 



C. Particle number dependence 



In Sect. IIV A[ we studied the Husimi distribution and 
the Wehrl-Husimi entropy for Yang-Mills quantum sys- 
tem by using TV = 1000 test particles. We note that 
the results of the test particle method we used to ob- 
tain Sji(t) depend on the number of test particles. The 
Husimi distribution pu{t](\, p) depends on the particle 
number N through the ansatz in (|41 [) . and so does the 
Wehrl-Husimi entropy S#(t). 

Our main goal in this section is to quantify the depen- 
dence of the saturated Wehrl-Husimi entropy on the test 
particle number N . We proceed with this study by the 
following method. First, we plot the energy histograms 
for several different numbers of test particles (we choose 
N = 1000, N = 3000 and N = 8000) in Fig. The 
distribution of the initial locations of the test particles 
in the phase space are generated according to the nor- 
mal distribution (f> defined in (|47|) . with the parameters 
given in (f57"l l60| . Figure [5] shows that the ranges of the 
test particle energies differ only slightly for N = 1000, 
N = 3000, and N = 8000. In other words, for the en- 
ergy distribution nxp (e) defined in (|64|). the center fi 
and width a are independent of N, notwithstanding small 
fluctuations. By fitting the energy histograms for various 



with fluctuations less than 0.5% and 5% respectively. We 
also define the normalized energy distribution of the test 
particles as 



n T p(e) 



n T p(e) 



Jo den TP (e) 



(79) 



Thus we conclude that the energy histograms for all 
choices of N correspond to a unique normalized energy 
distribution, riTp(e), which is unaltered by the time evo- 
lution and independent of N, provided that N is suffi- 
ciently large. 

Next, we compute the Wehrl-Husimi entropy Sn(t) for 
N = 3000 under the same set of initial parameters (|55F 
EH [60]) we used in Sec. lTVAl to calculate S H (t) for N = 
1000. We plot the Wehrl-Husimi entropy S H (t) for the 
two values of N in Fig. [5] We observe that the Wehrl- 
Husimi entropy S H (*) for N = 1000 and N = 3000 agree 
well for t < 2, but gradually diverges when t > 2. For 
both cases, the entropy begins to saturate at almost the 
same time, viz., t > 6.5. However, the saturation values 
are different: for N = 3000, Sjj(t) saturates to 8.1, while 
for N = 1000, S ff (t) saturates to 7.6. 

Based on the above results, we decided to analyze the 
saturation values of Sn(t) as a function of N. From Fig.[S] 
we conclude that the saturation is reached for t > 6.5, 
independent of how large iV is. We thus use <Sjj(10) as 
a proxy for the saturation value of Sn{t)- In Fig. [71 we 
plot Sh (10) for several different test particle numbers N 
and fit the curve by the function Sfi t (N), defined as: 



■S3 



(80) 



where S2, S3 and a are parameters determined by the fit. 
We obtain: 



s 2 = 8.73, s 3 = 76.4, a = 0.6115. 



(81) 



If our hypothesis is correct that £#(10) represents the 
saturation value of Su{t) for any N, this implies that the 
saturated value of Sn(t) approaches 8.73 for TV — » 00 for 
the initial conditions chosen for our numerical simulation. 



D. Canonical partition function and entropy 

We now consider the Yang-Mills Hamiltonian system 
in (]12fl for various classical ensembles. In the following 
numerical calculation, we use the same numerical param- 
eters as those specified in Sect. IIV Al We begin by ob- 
taining the canonical partition function and the canonical 
entropy for this system. We first determine the temper- 
ature of canonical ensemble of the Hamiltonian in (|12p 
that would be reached if the system would approach ther- 
mal equilibrium. This temperature can be obtained by 
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FIG. 6: The Wehrl-Husimi entropy S H (t) for TV = 1000 and 
N = 3000 respectively. In both cases, the test particles are 
generated at t — by the same set of initial parameters in 
(|55H57I I60|l . The Wehrl-Husimi entropies for both values of 
N agree well for t < 2, but gradually diverge for t > 2. Sii{t) 
for N = 3000 saturates to 8.1, while S H {t) for N = 1000 
saturates to 7.6. The saturation level is reached in both cases 
for t > 6.5. 



the following procedure. First of all, the total energy 
of the system is defined in (1201) and was evaluated nu- 
merically to be £ [HhPh] = 100.6 ± 0.5, as shown in 
Sect. I IV Al On the other hand, the canonical ensemble 
average of the Hamiltonian H 

(H)c = ~ [ rfr q , p H exp (-H/T) , (82) 

where T is the temperature and the partition function is 
defined as 



/oo 
dT^ p exp {-n/T) . 
-OO 



(83) 



We then fix 1(H) c to the total energy of the system 
£[UhPh]: 



£ [HhPh] = (H}c, 



(84) 



from which we determine the temperature T x of the 
equivalent canonical ensemble. 

When we try to evaluate (|84|) by substituting the 
Hamiltonian of the Yang-Mills system ([12")) into (|82p , we 
encounter a problem associated with the classical limit 
of the quantum system. The integrals over qi and q2 ex- 
hibit a logarithmic divergence owing to the special form 
of the potential V(qi, (72), which vanishes along the axes 
qy = and q2 — 0. A classical particle can therefore 
escape toward infinity in the hyperbolic channels along 



FIG. 7: £#(10) for several different test particle numbers N, 
indicated by the blue dots. We fit the curve by a fit function 
<Sfit(TV) defined in (1801) . The fit parameters are shown in the 
figure. 



the qi, qi axes [35j. In contrast, the escape of a quantum 
mechanical particle to infinity is forbidden by quantum 
fluctuations. The channels grow narrower as the particle 
moves away from the origin, and more and more energy 
is required to localize the particle in the direction or- 
thogonal to the channel. The uncertainty relation thus 
provides for effectively finite boundary conditions; as a 
result, the energy levels of the quantum system are dis- 
cretized [36j . 

Matinyan and Miiller [3^, [33] showed that this quan- 
tum effect could be accounted for in the semi-classical 
limit by adding a harmonic term to the Hamiltonian: 



n = 



1 

2^ 



(pi 



■Pi) 



1 



222 

9 9i9 2 



9 



2mT 



ill 



?2 2 ) 



(85) 



where the last term encodes the quantum correction. 
Thus, instead of inserting the classical Hamiltonian into 
(1121) . we apply the Hamiltonian with quantum correc- 
tions: 



H 



where 



1 

2m 



(pi 



pi) 



1 



9 2 <ll<ll + \muj 2 (q\ + ql) , (86) 



2„2 



2m 2 T' 



(87) 



The additional term arises from the commutator of the 
kinetic and potential energy in the semiclassical expan- 
sion of the partition function 
g = h = 1 we can now solve eq. 



37] 



After setting m 
for the equivalent 
temperature T x . We obtain T x w 67.1 and u ~ 0.0863. 
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Starting from the Hamiltonian (|86|) . we obtain the par- 
tition function for the canonical ensemble [38L w9| : 



Z(u,T) = 



mT 3 / 2 



2i:gh 2 



exp 



K 



/ 2 4 

, / m or 



Ag 2 T j u V 4.9 2 T 



where Kq(z) denotes the modified Bcssel function of the 
second kind. Since the free energy in the canonical en- 
semble theory is F = —T In Z and the entropy is given by 
Sc = ~dF/dT, the entropy of our system in the canon- 
ical ensemble is: 



Sc{u,T) = - 




+ 



The partition function Z diverges for uj = 0, and so does 
the canonical entropy Sc- Both divergences are cured by 
the quantum correction to the Hamiltonian (|86p . In view 
of the discussion above, we obtain the canonical entropy 



as S c (co,T x ) « 9.70. 



E. Microcanonical distribution and entropy 

In this section, we compare the late-time Husimi dis- 
tribution to the microcanonical distribution. Since the 
Yang-Mills quantum system is an isolated system, we an- 
ticipate that the Husimi distribution after equilibration 
would approach the microcanonical distribution. 

We obtain the appropriate microcanonical distribution 
by the following procedure. First, we construct the mi- 
crocanonical distribution in the test particle space by 



pmc (q, p) = 4 

a Jo 



de S [H H (q, P) - c] «tp (e) , (90) 



where tin (q, p) is defined in (|53|). e is defined in (|6T|) . 
nxp(e) is defined in (|79[) . and S is the normalization con- 
stant. We note that the initial energy distribution for our 
system is not strictly a delta function 8[H u(c[, p) — e], 
because we generated the test particle positions in phase 
space randomly according to the distribution defined 
in eq. (jT7) . Therefore, Pmc (q, p) must be defined as 
<5[?iir(q,p) — e] folded with the energy distribution of 
test particles shown in (|90|). According to (|M|) . the en- 
ergy is conserved for each test particle individually, and 
thus nxp (e) remains unchanged as time evolves. Using 
([oil). ([791 and ([90]). we easily obtain 



1 



Pmc (q, P) = ^7 e x P 



1 - 2 
7^2 P) 



(91) 



where /i and cr are input from (|78|l . S' is the redefined 
normalization constant and T-Lh (q, p) is obtained from 
(|62|) . In the test particle space, pmc is normalized as: 



drq i pp M c(q, p) = 1. 



(92) 
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FIG. 8: Energy histogram of test functions for pmc (q>p), 
which is defined in (|91[) . The test functions are generated by 
Metropolis-Hastings algorithm, and the total number of test 
functions is M = 8 x 10 4 . e denotes the test-particle energy, 
which is defined in (|61[) . and the labels on the vertical axis 
denote test particle numbers. A normal distribution riMc(e) 
is used to fit this histogram. A, fiMC and <tmc are the fit 
parameters for rc.Mc(e), which are defined in (|64[) . The values 
for the fit parameters are shown in the plot. 
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FIG. 9: u-histogram of test functions for pivic (q, p), which 
is defined in (|91[) . The test functions are generated by 
Metropolis-Hasting Algorithm, and the total number of test 
functions is M — 8 x 10 4 . u is defined as u — qi<]2, the labels 
on the vertical axes denote test-particle numbers. 
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To obtain the microcanonical distribution in the phase 
space PMc(q ; p)i we convolute pmc with test particle 
Gaussian K, which yields: 

Pmc (q, P) 

/OO 
dTq,pPMC (q, p) K (q - q, p - p) , (93) 
- OO 

where pmc is defined in (|9"Tj) and K is defined in (j4l)|). 
The microcanonical entropy is then obtained as: 



Smc — — 



<ff 1 q, P /'Mc(q,p)lnjOMc(q,p)- (94) 



Before we proceed, we briefly comment on the reason 
why pmc (q>p) should be constructed by (|93l) . In statis- 
tical physics, the microcanonical distribution of an iso- 
lated system of energy E is conventionally obtained by 
Pmc = S(T~L — E)/ft, where ft is the total number of 
microstates that satisfies the constraint H = E. If we 
substitute this conventional definition of pmc into (pHf. 
it is straightforward to show that Smc is not well defined. 
However, if one approximates 5{T-L — E) by an Gaussian 
distribution centered on E with a finite width <j g , Smc 
becomes well defined and is a function of both, E and a g . 
Therefore, pmc( x 7P) in (1531 is defined in a way that en- 
codes the coarse grained energy of the system, the width 
of energy distribution and the widths for the test parti- 
cle Gaussians, all of which must be equivalent to those 
specified in our choice of the initial Husimi distribution 
p ff (0;x,p) in Sect. HUH and W~M 

Owing to the complexity of (|91[) and the multidimen- 
sional integrals (|93l) and (|M|) . we adopt an alternative 
approach to evaluate pmc (q, p) , instead of directly eval- 
uating eq. Our approach is briefly described as 
follows. Since pmc (q 5 p) in (ED is a non-negative func- 
tion and normalized by (|92p . we generate a sufficiently 
large number of test functions in (q, p)-space according 
to the distribution pmc (qp)- Thus pmc (qp) can be 
represented as a sum of these test functions: 



1 M 

PMc(q, p) = ]jj £ [«J(q - q")<*(p - p s )} 



(95) 



where (q s , p s ) denotes the locations of the test functions, 
and M is the total number of test functions. We gener- 
ate (q s ,p s ) by the Metropolis-Hastings algorithm using 
5 x 10 6 iterations. After excluding the first 10 5 itera- 
tions, we randomly select, for instance, M = 8 X 10 4 
points (q s ,p s ) from the remaining 4.9 x 10 6 iterations. 
In view of the shapes of the position and momentum pro- 
jections of pmc (q, p) , we make the following change of 
coordinates: u = qiq 2 and v — tan -1 ^)- To ensure that 
the locations of the test functions are ergodic in (q, p)- 
space, we impose periodic boundary conditions to the 
random walks in the Metropolis-Hastings algorithm. For 
instance, when setting p = 100.6 and a = 8 in (|iJTj) . we 
can map the entire domain in each dimension periodically 



to the region: \u\ < 16, \v\ < (tt/2 - 10~ 5 ), |pi| < 16.5 
and \pz\ < 16.5. In this case, the acceptance rate is about 
22%. 

To verify the validity of the resulting microcanonical 
distribution, we plot the energy histogram of the test 
functions and compare it to the energy histogram of the 
test particles used to represent the Husimi distribution. 
In Fig. [51 we plot the energy of for the test functions 
for the microcanonical distribution. According to (|61l) . 
e s = Hh (q s , P s ) denotes the energy for the test function 
s, for s = 1,...,M. We fit the energy histogram for the 
test functions for pmc (q, P) by the normal distribution 



riMC (e) = A exp 



2a 



1 f \» 
2 — (.£ — MmcJ 

MC 



(96) 



The values of the fit parameters A, pmc and <tmc arc 
listed in Fig. [8] for M = 8 x 10 4 . We obtain: 



Pmc = 101.1, ctmc = 7.975. 



(97) 



We define the normalized energy distribution for test 
functions as 



«Mc(e) 



fiMc(e) 



J den M c(e) 



(98) 



Comparing ([75)1 to (|57)) . we obtain pmc ~ A* an d ^mc ~ 
a, with the errors less than 0.5%. Therefore, we conclude 
that nMc(f) in ([98]) is practically identical to nTp(e) in 
([79"lh with the errors of less than 0.5%. Furthermore, 
in Fig. [9] we plot the w-histogram of the test functions 
for pmc (q, p), where u = q\q 2 . Figure |9] shows that 
the distribution of test functions is symmetric in the u 
coordinate. 

Substituting (f9"5j) to ([SBJ, we obtain: 



Pmc 



1 M 



8=1 



where K is defined in (j46]l and we choose 7^- = 3/2 in 
(ISUl) . Clearly, pmc is normalized by: 



/OO 
dr q ,pPMc(q,p) = 1. 
-00 



(100) 



We visualize pmc (q, p) in (|99|) by projecting on the 
(#1)92) and {p\,p 2 ) subspaces, respectively: 

/OO 
dpidp 2 p M c(gi,<72,Pi,P2),(101) 
-00 

/OO 
dqidq 2 Pmc (?i, 92,Pi,P2)-(102) 
-OO 

We plot F q ylc (q 1 ,q 2 ) and F^[ PllP2 ) in Fig. HOI When 
we compare Fig. [3] to Fig. I10[ we find that F q and F p 
at time t = 10 are very similar in shape to F q ylc and 
F? 4C , respectively. Contour lines of both F q (t = 10) and 
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FIG. 11: Comparison of G(t;p) at t = 10 and Gmc(p)- 
We define G(t;p) and Gmc(p) in (1103|) and (1104|) respec- 
tively. G(10;p) is obtained from the momentum projec- 
tion of ph(W; q, p) composed of N = 10 4 test particles, 
while Gmc(p) is obtained from the momentum projection of 
PMc(q, p) composed of M — 2 x 10 4 test functions. 

F^ 1 follow equipotential curves, while the contour lines 
of both F p (t = 10) and are shaped as concentric 

circles. 

To quantify the similarities between pn(t] Q, p) at late 
times and pmc (Qi p) i we compare their momentum pro- 
jections. By switching to polar coordinates p\ = pcos9 



20 r j 




2 -% -15 -10 -5 5 10 15 20 



microcanonical distribution function (a) J 5 ^ 10 (91,92) and (b) 
are generated by Metropolis-Hastings algorithm, and the total 



and P2 — psmO, we define the following two projections: 

G(t;p) = / dO F p (t; p cos 6, p sin 9) , (103) 
Jo 

G M c(p) = / d0 F p MC (p cos 0, p sin 9), (104) 
Jo 

where F p and F^ c are defined in ([66]) and (|102[) respec- 
tively. In Fig. [TTJ we plot G(10;p) and Gmc(p) m com- 
parison. G(10; p) is obtained from the momentum projec- 
tion of pff(10;q, p) composed of N — 10 4 test particles, 
and Gmc(p) is obtained from the momentum projection 
of PMc(q,p) composed of M = 2 X 10 4 test functions. 
The figure shows that G(10;p) and Gmc(p) have similar 
values for all p, and the largest deviation occurs at low 
p. G(10;p) and GmcCp) a ^ l° w P receive contributions 
from the test functions located at the narrow "channels" 
along the coordinate axes in the position projections of 
Ph and pmc, respectively. Since the number of test func- 
tions, N and M, are finite, one expects larger fluctua- 
tions of the contributions from these narrow "channels" , 
which explains the observed deviation at small p. Over- 
all, the close similarity between G(10;p) and Gmc(p) sug- 
gests that ;0#(i;q, p) asymptotically approaches the mi- 
crocanonical density distribution pmc (q, p) ■ 

Finally, we obtain the microcanonical entropy Smc by 
substituting (|M| into (|M|). We evaluated Smc with the 
help of Simpson's rule and by applying the same inte- 
gration techniques as those described in Sect. IIV Al We 
verified the numerical precision of our approach by evalu- 
ating the normalization for pmc(q,p) for various choices 
of M and found that the numerical result coincides with 
(IIOOP within errors of less than 0.6%. In addition to the 
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FIG. 12: The microcanonical entropy Smc as a function of 
M, indicated by the blue dots. Smc is denned in (I94|) . M 
denotes the total number of test functions, as revealed in (|95[) 
and (J99J) . We set y, = 100.6 and cr = 8 in (J9TJ. Besides, we 
fit the curve by a fit function Sm(M) defined in (|105[l . The 
fit parameters are shown in the figure. 

errors associated with the use of Simpson's rule, Smc pos- 
sesses an additional error, typically less than 0.5%, which 
arises from the Monte-Carlo calculation of /9Mc(q, p) in 
([9"5"j). In Fig. Q21 we plot Smc for several different test 
function numbers M. We fit the data by the function 

S fit (M) = S4 -^. (105) 

The parameters determined by the fit are: 

s 4 = 8.788, s 5 = 1258, c = 0.9517. (106) 

We thus conclude that Smc ~ 8.79 is the microcanonical 
entropy for our chosen initial conditions. 

In Sect. HVA[ we obtained the value Sn(t = 10) -> 
8.73 in the limit N — > oo for the initial conditions cho- 
sen for our numerical simulation. Under the same initial 
conditions, we found Smc — > 8.79 when M — > oo. We 
conclude that the saturation value of the Wehrl-Husimi 
entropy coincides with the microcanonical entropy within 
errors, estimated at 1%. Apart from numerical errors, the 
difference between the two entropy values may also be ac- 
counted for by the fact that at t — 10 the system may 
not yet be completely equilibrated. Since Smc < Sc, 
we also conclude that the Yang-Mills quantum system is 
equilibrated microcanonically but not thermalized. The 
system does not have enough degrees of freedom to render 
the microcanonical and the canonical ensemble approxi- 
mately identical. 

In the above calculation, we studied the micro- 
canonical distribution Smc fo r the Yang-Mills quan- 
tum mechanics model at the coarse grained energy \i = 



FIG. 13: The microcanonical entropy Smc as a function of 
M for the coarse grained energies fi = 50.6, 100.6, and 200.6. 
The corresponding widths a, defined in (|9f [) . for these energies 
are a = 5.8, 8.0, and 11.5. We fitted these points by the 
function Sfm(M) defined in (|f05|) . and use the fit parameters 
to determine the asymptotic values of Smc for M — > oo, which 
are Smc = 7.88, 8.77, and 9.54 (from bottom to top). 

£[HhPh] ~ 100.6. We now briefly comment on how Smc 
depends on the coarse grained energy of the system. In 
Appendix[C] we show that while the Yang- Mills Hamilto- 
nian % possesses a scale invariance, the scale invariance of 
%h is partially broken when we demand that the smear- 
ing function in ([5]) should retains its minimal uncertainty. 
The reason is that, for any coarse grained average energy 
H, the relation £7? = ti 2 /4 constrains our ability to rescale 
£ and r\ in (fl"8|) . Alternatively, we observe that the ad- 
ditional terms in the expression for Hh (q, p) break the 
scaling symmetry of the original Yang-Mills Hamiltonian. 

Despite the fact that the scaling properties of Hh are 
partially broken, we can examine numerically how Smc 
changes when \x scales as /i — > A*//, where X s is the scal- 
ing parameter. In analogy to (|C4|) . we parametrize the 
change in the microcanonical entropy as 

Smc (a 1 ) -> Smc(^) + rln\ s , (107) 

where r is a constant to be determined numerically. In 
order to find the value of r, we calculated Smc by nu- 
merically evaluating for various choices of /x in (|9ip . 
In Fig. [T31 we show Smc a s a function of M for /1 = 50.6, 
fi = 100.6, and fj, = 200.6, respectively. The corre- 
sponding widths a, defined in (f9"Tj) , for these energies are 
a = 5.8, 8.0, and 11.5, respectively. In Fig.QSl we fitted 
these curves by S&t(M) defined in (|105[) . The fit param- 
eters again determine the asymptotic values of Smc for 
M -> 00. The results are S M c = 7.88, 8.77, and 9.54, 
respectively. From these results we can deduce the value 
r = 5.0 ±0.2. 
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In Appendix [C] we show that the scale invariant Yang- 
Mills Hamiltonian % implies the value r' = 6, where r' 
is defined in (|C4j) . The difference between r and r' is 
attributed to the following reason: Since we demand the 
Gaussian smearing function in ([5]) retains its minimal 
uncertainty encoded in the relation £77 = h 2 /4, we are 
breaking the scaling symmetry of the Husimi Hamilto- 
nian "Hh, as discussed in Appendix [Cl This argument 
suggests that Smc (a«) changes less strongly under a scale 
transformation than naively expected. Comparing the 
numerical value for r with the analytical value for r' , we 
indeed obtain r < r' , which confirms our expectation. 

V. CONCLUSIONS 

We have developed a general method to solve the 
Husimi equation of motion for two-dimensional quan- 
tum mechanical systems. We proposed a new method 
for obtaining the coarse grained Hamiltonian whose ex- 
pectation value serves as a constant of motion for the 
time evolution of Husimi distribution. We solved the 
Husimi equation of motion by the Gaussian test particle 
method, using fixed-width Gaussians. Having obtained 
the Husimi distribution, we evaluated the Wehrl-Husimi 
entropy as a function of time for the Yang-Mills quantum 
system. 

By comparing the Wehrl-Husimi entropy Sn(t) ob- 
tained from different particle numbers, N = 1000 and 
N = 3000, we found that the values of Sji(t) agree 
for t < 2, and saturation is reached in both cases af- 
ter t > 6.5. However, Sji(t) for N = 3000 saturates to 
a higher value than for N = 1000. This result suggests 
that for a larger number of test particles the Husimi dis- 
tribution is more evenly distributed in the phase space, 
and thus a larger value of N results in a higher satu- 



ration value of the Wehrl-Husimi entropy. By evaluat- 
ing Sff(10) for various different A's, we concluded that 
Sh (10) — > 8.73 for N — > 00 for our chosen initial condi- 
tions. 

In order to address the question of equilibration, we 
studied the Yang-Mills Hamiltonian system in the canon- 
ical and microcanonical ensembles. The canonical en- 
tropy for the system is Sc ~ 9.70. We obtained the 
microcanonical distribution by generating M test func- 
tions. We observed that the saturated Husimi dis- 
tribution closely resembles the microcanonical distribu- 
tion. Moreover, we obtained the microcanonical entropy 
Smc —> 8.79 as M — > 00 for the same choice of initial con- 
ditions. Therefore, comparing the saturation value of the 
Wehrl-Husimi entropy to the microcanonical and canon- 
ical entropies, we conclude that (SH)max ~ Smc < Sc- 
This implies that, at late times, the Yang-Mills quantum 
system is microcanonically equilibrated but not thermal- 
ized. 

It is straightforward to generalize the method intro- 
duced here to solve the Husimi equation of motion in 
three or more dimensions. However, for higher dimen- 
sions, the evaluation of the Wehrl-Husimi entropy be- 
comes even more challenging owing to the increasing 
numbers of integrals. A new method will then be needed 
for the reliable evaluation of entropy. 
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Appendix A: Equations of motion for the test particles 



In Sect. Hill we obtained the equations of motion for the ten variables, q[, q%, p\, p\, ct igi , c* 2?2 , ci ipi , c* 2 , c l qipi 
and c l q2P2 , where i labels the test particle. In (|2"51 lUT]) . we listed the equations obtained from the first moments I qi , 
I Q2 , I Pl and I P2 of eq. (|13p. The equations obtained from the second moments I q 2, I q 2, I p 2, 7 p |, I qiPl , and I q2P2 of 
(fT3")) are listed below: 



- O*) ( c U (*)) 2 - O*) KpM 2 + -4 lP1 (*) A iw = 0, 



w q2P2 (ty q2P2 (t)c; 2P2 (t) 



.(*) -<W*) 



c q2P2 (t)AUt) = 0, 



(Al) 
(A2) 
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24 P1 (*)0*)4*(*)- 


O*) (4» (*))'-&»(*) («U(*)) a 


- 


_ 




/ 4 2p2 (t) aA 












dqldql 


qi(t). 






















9g| 


*(«) ^ 2 y 


dq\dq\ 


qi(t)_ 



4 m (i)Aj(i) = 0, 



4 2P2 (i)A* 2 (*)=0, 



4« (*)4x PI (*)4 P x w + 4 P1 (*)4* (*)4» (*) - 4*, (*) (4« (*)4 K (*) + (4« (*)) : 

a 2 



2ma to I A\(t) I \ A\(t) 2 J dq\ 



q*(t) 



2 I Ai(t) 



2 A|(i) 2 / 9g 2 9g| 



q«(t). 



(Ai(f)) =0, 



'92 92 



(*)4»(*)4,»(') + O^^UW - (4*(')4«»(*) + (4 2P3 (*)) ; 



2to« to \ A^i) 

i few 



A l 2 W 2 a J 9g 2 



q*(t) 





d 4 V 




/ ^ A|(t) 2 ) 


dqldql 


q«(t)_ 



2 ^ AKO 

where i = 1, 2, 2V, and Aj (t), A 2 (t) and q l (t) are denned in ([2"5]l. (|2l)| and ([52]). respectively. 



(A^)) 2 =0, 



(A3) 



(A4) 



(A5) 



(A6) 



Appendix B: Husimi equation of motion in one dimension 

The Husimi equation of motion for one-dimensional quantum systems was derived in [l9| . For the potential energy 
V(q) being a C°°-differentiable function of g, the Husimi equation of motion in one dimension is: 



dpH 
dt 



1 



h 2 d\ d PH 

2a dp J dq 



A-l 



d x+f *V{x) d x d^- 2K 



2 a+ m -i x\K\(fi- 2«)! dg A +^ dp x dq^- 2K 



PH 



(Bl) 



where A, p and « are summed over all non- negative 
integers subject to the constraints that A is odd and 
p - 2k > 0. 

We discuss the energy conservation for the one- 
dimensional Hamiltonian. As a specific example, we start 
from the following one-dimensional Hamiltonian in the 
Wigner representation. 

*<*■'> = £-f' 3 + ii» 4 ' (B2 > 

where A and £ are positive-valued parameters. We de- 



rive the corresponding one-dimensional coarse grained 
Hamiltonian as follows. The Husimi distribution for a 
one-dimensional quantum system can be obtained from 
the Wigner distribution by: 

p H (t;q,p) = — / dq'dp' e -(<?'-9) 2 /«-a(p'-P) 2 M 2 

xW(t;q',p'). (B3) 

Starting from (|B3[) and proceeding like in Sect. Ill C[ we 
obtain an expression similar to that of (|17p . which reads 



Hh (q,p) 



(27T) 



dx'dp' %{q' ,p') I dudvexp 



— u 2 + - — v 2 — iuiq 1 — q) — iv(p — p) 
4 4a 



(B4) 
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FIG. 14: Solutions of the Husimi equation of motion in one dimension. The Hamiltonian is defined in (|B2[I of Appendix iBl 
The parameters as chosen as re = 1 and £ = 1. The initial conditions are discussed in Appendix [B] Panels (a) and (b) show 
pH{t;x,p) for a single test particle, at time (a) t = and (b) t — 2. Panels (c) and (d) show pH(t;q,p) or the many test 
particles, at times (c) t = and (d) t=2. It is obvious that for t > this single particle ansatz is insufficient to represent the 
solution. 



Here u and v are Fourier conjugate variables to q and 
p respectively. Similar to the calculation in Sect. Ill CI 
we set £ = —a/4 and r\ = —h 2 /(Aa). We evaluate the 
integrals in (|B4j) in the analytic region where £ > and 
i] > 0, and then substitute £ = —a/4 and ?y = — h 2 /{Aa) 
into the resulting analytical expression. In this manner, 



we obtain the coarse grained Hamiltonian: 

~^ + h a ^ +8K) - (B5) 



Proceeding similarly as in Sect. Ill CI we use eqs. (|B1[) 
and (|B5p to prove that £ [HhPh] is a constant of mo- 
tion for the Husimi equation of motion in one dimension. 
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Thus S [HhPh] should be identified as the total energy 
corresponding to the Hamiltonian (|F32|) . 

Next, we solve the Husimi equation of motion (|Blj) by 
using the test-particle method described in Sect. IIIII We 
begin by writing the Husimi distribution as: 



Applying the five moments I q , I p , I q 2 , I p 2 and I qp to the 
Husimi equation of motion (|B1|) . we obtain five equations 
of motions for each test particle i for the five variables 
representing the location in phase space and width of 
each test particle. 



p H (t;q,p) = — V Az (t) exp 



i=l 



x exp 



1 

2 



c PP (t) ( P -i?(t)) 



X exp [-c qp (t) (q - q\t)) (p - p\t))} , (B6) 
where i = 1 , N, and we define 



= t qq {ty pp {t) - (c qp (t)y 



(B7) 



The moment of a function f(t;q,p) with respect to a 
weight function w(q,p) is defined as: 



I W [f] 



dq dp 
2-Kh 



w{q,p)f(t;q,p)} 



(B8) 



These equations are: 



dV 



. A* (t) 



- -p\t) 
m 



a 
2 



93 y 



9g 3 



0, (B9) 



= 0, (BIO) 



9*W 



2c* p (i)c* p (i)4 p (t) - c qq (t) (c pp (t)) 2 - ^(t) 



-c^(i)A J (i)=0, 



(Bll) 



^qp C qp C qq ^qq W ( C <?p W) ^pp(0 ( C ijij(^)) 









<9 4 ^ 




[ 9g 2 


+ 

g*(t) 


^A*(i) 2 J 


dg 4 


r(t)_ 



4,(t)A^) = o, 



(B12) 



qq\"J~qp\ 

h 2 i_ 

m 



2ma 



A*(t) 



A*(^ 
\ 2 



--a 



a 
2" 



d 4 v 



<9 9 4 



g-(t) 



A'(t) 
(A l (t)) 



9<y 2 



(B13) 



where i = 1, iV. By solving (|F39|) - (|B13I) simultaneously 
for j = l,...,N, 
functions of time. 



for i = 1,...,N, we obtain q\ p 1 , c z xx , c pp and c l xp as 



Finally, we solve these 57V equations of motions for the 
Hamiltonian system in (|B2|) . with k = £ = 1. For choos- 
ing the initial conditions, we adopt the method similar 
to that introduced in Sect. IIII 51 Here, we briefly outline 
the ideas without showing the details. We choose the ini- 
tial conditions setting the initial Husimi distribution to 



be: 

f°° dq'dp' , 
p H {0;q,P) = / -——-K{q-q,p-p) 
J-00 27rft 

xcf>(q',p'), (B14) 

where K and <j> are defined in Sect. IIII Bl We express pH, 
K and <f> in the forms of (H31 . (|46p and (|4"T)) respectively, 
with the redefined variables x — (q-,p) an d x' = W ,p') 
and the redefined indices a — 1,2 for x a : x' a > Ph> P%-> 
1hi Ik an( i it- -^y convolution theorem, we obtain that: 

— = — + — , (B15) 

l a H lK 1% 
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for a = 1,2. At t = 0, we choose 7^ = 1. In the 
many-particle ansatz, we choose N — 1000, 7^- = 3/2 
and 7^ = 3. And we choose p H — {A — 0. In the sin- 
gle particle ansatz, pn remains a single Gaussian for all 
times, and thus we choose jfj — 1 and p a H = 0. We plot 
Pii{t;q,p) both for the single-particle and many-particle 
ansatz in Fig.[TJ] We discuss the meaning of these results 
in Sect. IITTB1 

Appendix C: Effects of coarse graining on the scale 
invariance of the Yang-Mills Hamiltonian 

In this section, we discuss the effects of coarse grain- 
ing on the scale invariance of the Yang-Mills Hamilto- 
nian. We begin by constructing an alternative micro- 
canonical distribution p' M q in terms of the conventional 
Hamiltonian H in (fl~2f and the conventional energy E, 
and we obtain the scaling of the microcanonical entropy 
»!?mc with respect to that of E. Furthermore, we show 
that, while % is scale invariant, the scale invariance of 
the coarse grained Hamiltonian Hh is partially broken, 
due to the requirement that the smearing Gaussian func- 
tion in the Husimi transformation ([5]) should retain its 
minimal quantum mechanical uncertainty. 

For the conventional Hamiltonian in (|12p , we construct 
an alternative microcanonical distribution p' M q as: 

^c = ^ex P (-^J. (CI) 

As discussed in Sect. lIVEl approximating 6(H — E) by a 
Gaussian distribution is a way to construct a microcanon- 
ical distribution that leads to a well-defined entropy. De- 
fine A s as a scaling parameter. As the position and mo- 
mentum scales as q — > A s q and p —> A^p respectively, 
it is straightforward to show that T-L — > and thus 

E — > XgE. The normalization condition: 

/rfr q , p p^ c (q,p) = l (C2) 



must be scale invariant. Owing to the scaling r q p — > 
Afr qiP we obtain £1 — > A^fi and a g — > \ A s <J g . The micro- 
canonical canonical entropy S^c * s defined as: 

S'mc = -J dT^ p p' MC {q,p)\np' MC (q,p), (C3) 

where p' MC is given in (|C1|) . The scaling of S^ c follows 
from the scaling of H and E: 

S^ c (E)^S^ c (E) + r'\n\ s , (C4) 

where r' = 6. 

The coarse grained Hamiltonian Hh(<1, p) given in (fH))) 
is obtained from H(q, p) by the transformation (JTS)). We 
now examine how Hh (q, p) scales when the positions and 
momenta scale as q — » A s q and p — > A^p, respectively. 
The uncertainty relation of a quantum state reads: 

AqiApj^-Sij, (C5) 

where i,j = 1,2. We note the difference by a factor of 
2 between (|C5[) and (|35|) . which was pointed out in |25| . 
From (fTB)) and (|C5|) . it is straightforward to show that, 
when q — > A s q and p — > A^p, Hh will scale as Hh — > 
\ a s T-Lh only if the smearing parameters £ and rj scale as 
£ — > and 77 — > X^r), respectively. In addition, the 
constraint A s > 1 is imposed by the uncertainty relation 

The Husimi distribution is defined as a minimally 
smeared Wigner function, as can be seen from ([5]). For 
the smearing Gaussian with minimal uncertainty, we 
have AqjApj = h/2 for j = 1,2, and thus £77 = h 2 /4. 
Therefore, we do not have the flexibility to scale the pa- 
rameters £ and 77 in the required way, if we demand that 
the smearing Gaussian in (J5J) should retains its minimal 
uncertainty. As a consequence, the scaling symmetry of 
Hh is partially broken. We compare the different scaling 
behavior of S^ LC (E) and SmcC/- 4 ) & t the end of Sect. lIVEl 
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